import IPython.core.display as di
# This line will hide code by default when the notebook is exported as HTML
di.display_html('<script>jQuery(function() {if (jQuery("body.notebook_app").length == 0) { jQuery(".input_area").toggle(); jQuery(".prompt").toggle();}});</script>', raw=True)
# This line will add a button to toggle visibility of code blocks, for use with the HTML export version
di.display_html('''<button onclick="jQuery('.input_area').toggle(); jQuery('.prompt').toggle();">Toggle code</button>''', raw=True)
#!pip install folium
#%load_ext autoreload
#%autoreload 2
import matplotlib
import matplotlib.pyplot as plt
import pandas as pd
import cufflinks as cf
import plotly.offline as py
import plotly.graph_objs as go
from IPython.display import HTML
import numpy as np
import seaborn as sns
#! pip install hex_bin_maker
from hex_bin_maker import HexBinMaker
#OI IP
! pip install workalendar
import OFL.base.oi_colors as oic
from proteus_ds import RawObjMsmt
from OFL.base.utils import geom as oi_geom
from OFL.base.dataframe import shapes as oi_pandas
import os
!pip install folium
import folium
from folium import FeatureGroup, LayerControl, Map, Marker
%matplotlib inline
py.init_notebook_mode()
cf.go_offline()
cf.set_config_file(theme='white')
Patterns of Life: Marawi City, Philippines¶
Automated Monitoring of Car Counts, Landcover Change, and Geolocation Observation¶
Orbital Insight leverages proprietary imagery ingestion and computer vision techniques to count cars and trucks and analyze land use/land cover (LULC) and change detection at scale. In addition, we have begun to incorporate geolocation data into our analysis to create a more accurate portrayal of pattern of life, anomalies, etc.
This analysis helps analysts better understand what is happening on the ground in areas that are otherwise non-permissible. Orbital Insight Public Sector Solutions Engineers leveraged proprietary algorithms and data sources to understand pattern of life trends before, during, and after the Battle of Marawi.
This project shows some general tools and workflows that an analyst could easily deploy to produce information on a conflict area. First, car counts in the most heavily affected area of the city will be shown. Next, data from our landuse algorithm will be used to create visuals that show changes in building and road land cover classes. Lastly, sample geolocation data for the Marawi area will be shown against that of a nearby city where fighting did not take place.
Car Counts¶
The use of Orbital Insight’s car counting algorithm was first used to establish a pattern of life (PoL) in the most affected area of Marawi. Resultant trends identified through our object detection and data science/analytics can enable automated analytics and trend monitoring for:
- Major events on the ground.
- Removal of cars and people in an otherwise populated region.
aoi_dict = [{'city':{'Marawi_small':'Marawi'
}
}
]
raw_counts_city = RawObjMsmt(proteus_project_ids=['marawicarsearlie-lgSHIj'],
aoi_hierarchy=aoi_dict,
aggregation_level='city',
max_pct_aoi_unobserved=0.5,
min_obs_pct=0.9,
max_cloud_score=0.5,
token='be8929ebce4e1425cdd9494323c14e01e9a139d0')
raw_df = raw_counts_city.raw_df
plt_df = raw_df.groupby(['city'])['impute_status'].value_counts().reset_index([0])
plt_df['count'] = plt_df['impute_status']
plt_df['impute_status'] = plt_df.index
plt_df = plt_df.reset_index(drop=True)
plt_df = plt_df[plt_df['impute_status']!='not enough observed tiles']
plt_df = plt_df.pivot(index='impute_status',columns='city',values='count').T
Available Imagery and Processing Metrics¶
DigitalGlobe high-resolution imagery was the sole provider utilized for car detection. This was determined through the high performance and relevance of the DG Car Counter to this specific analysis. Airbus imagery could also have been explored to supplement some of the more sparsely collected AOIs.
Once imagery is ingested into our system, and tiling and pre-processing are completed, algorithms can be processed relatively quickly. This varies by project depending on how many square kilometers are undergoing analysis.
For the area most affected by the fighting, approximately 148 scenes were analyzed where 35,00 cars were counted for. This process took approximately 15 minutes of human analyst time and approximately 24 hours of compute time to order, ingest, and analyze all imagery. The larger the AOI, the longer it takes for all required imagery to be ingested. Also, note that due to the geography of the region, many scenes were unusable due to cloud cover.
from IPython.display import Image
from IPython.core.display import HTML
Image(url = "https://dl.dropboxusercontent.com/s/swu7tjav94jyrgi/marawi_small_hexbin_cars.png?dl=0?raw=1", width = 800)
Marawi Car Counts¶
Summary and Workflow:¶
The Battle of Marawi was a five month long conflict that started 23 May 2017 and is considered the longest urban battle in modern Philippine history. The AOI examined is the southeastern part of Marawi, which is separated from the "safe zone" to its northwest by the Agus River, and bordered to its south/southwest by Lake Lanao.
This AOI is where the heaviest fighting and most destruction occured, and therefore possible changes in car counts are expected.
Analytic Highlights:¶
Based on the time series of this AOI, we note that raw car counts immediately decline once the battle starts on 23 May 2017. These trends are consistent with the expectation that commercial and passenger traffic through Marawi are serverly disrupted due to the battle. The severe drop in automated car counts as soon as the battle begins would warrant a more in-depth analysis. A user can set up alerts for car count increases and decreases, and be notified if a specific threshold is met. Upon receiving a "decline in car counts" notification, an analyst would be prompted to examine the area more closely.
Main takeaways include:
- Notable drops in raw and rolling mean values during the initial start of the battle in late May 2017.
- As of July 2018, car counts have not completely recovered to their pre-battle numbers, possibly indicating that the city is still rebuilding, lacking infrastructure, and not suitable for repopulation. However, there is some hint of car numbers starting to increase after the battle possibly indicating cleanup and reconstruction. These numbers were taken from the available DigitalGlobe imagery dating back to the beginning of 2009.
*Due to the relatively large size of the AOI, not all scenes provide 100% coverage. A threshold for the percentage of AOI unobserved can be set manually. This allows for the customization between using more scenes and car count observations, even if the scene only covers a small part of the AOI, versus less scenes that may have a more complete coverage profile. Therefore, it is likely that the very small counts shown are the result of limited coverage scenes, or scenes that are only covering areas that normally have low car counts (a field versus a parking lot). In this case, the rolling mean of car counts can be used as a better indicator of PoL. Average car counts per square kilometer was there also used as another way of verifying that car counts did in fact decrease once the battle began.
df = pd.read_csv("/orbital/use_cases/marawi/count_per_area.csv")
data = [
go.Scatter(
x=df['imaging_ts'],
y=df['count_per_sqkm'].rolling(10).mean(),
name='Rolling Mean',
mode = 'lines+markers',
line = dict(color = (oic.OI_GREEN_HEX),
width = 2),
marker = dict(color = (oic.OI_GREEN_HEX),
size = 5)),
go.Scatter(
x=df['imaging_ts'],
y=df['count_per_sqkm'],
name='Marawi',
mode = 'lines+markers',
line = dict(color = (oic.OI_DARK_BLUE_HEX),
width = 2),
marker = dict(color = (oic.OI_DARK_BLUE_HEX),
size = 5)
),
go.Scatter(
x=['2017-08-10'],
y=[0,0,0,0,0,0,0,0,0,0],
mode='markers',
text=[ 'Battle of Marawi'
],
marker = dict(
color = '#d3d3d3',
size = 5,
opacity = 0.0,
line = dict(
color = 'black',
width = 2
)
),
showlegend = False
),
]
layout = go.Layout(
title = 'Average Car Counts Per Square Kilometer',
xaxis = dict(
title='Date'),
yaxis = dict(
title='Cars Per SqKm'),
# to highlight the timestamp we use shapes and create a rectangular
shapes = [
# 1st highlight during 23 May 2017 - 17 Oct 2017
{
'type': 'rect',
'layer':'below',
# x-reference is assigned to the x-values
'xref': 'x',
# y-reference is assigned to the plot paper [0,1]
'yref': 'paper',
'x0': '2017-5-23',
'y0': 0,
'x1': '2017-10-17',
'y1': 1,
'fillcolor':oic.OI_LIGHT_BLUE_HEX ,
'opacity': .5,
'line': {
'width': 0,
}
}
]
)
fig = go.Figure(data=data, layout=layout)
py.iplot(fig,filename='cufflinks/line')
plt_df_ = raw_df[raw_df['pct_aoi_unobserved']<0.75]
plt_df_ = plt_df_[plt_df_['city'].isin(['Marawi'])]
mean_L1 = plt_df_.groupby(['rounded_local_time','city'])['mean_count_L1'].sum()
plt_df = plt_df_.groupby(['rounded_local_time','city'])['obj_count'].sum().reset_index()
plt_df['rolling_mean'] = plt_df['obj_count'].rolling(window=10).mean()
plt_df.index = plt_df['rounded_local_time']
plt_df['mean_count_L1'] = mean_L1.values
plt_df = plt_df['2009-1-23':]
data = [
go.Scatter(
x=plt_df['rounded_local_time'][plt_df['city']=='Marawi'],
y=plt_df['obj_count'][plt_df['city']=='Marawi'],
name='Marawi',
mode = 'lines+markers',
line = dict(color = (oic.OI_DARK_BLUE_HEX),
width = 2),
marker = dict(color = (oic.OI_DARK_BLUE_HEX),
size = 5)
),
go.Scatter(
x=plt_df['rounded_local_time'][plt_df['city']=='Marawi'],
y=plt_df['rolling_mean'][plt_df['city']=='Marawi'],
name='Rolling Mean',
mode = 'lines',
marker = dict(color=(oic.OI_GREEN_HEX),
size=5),
line = dict(color = ('oic.OI_GREEN_HEX'),
width = 2)
),
go.Scatter(
x=['2017-08-10'],
y=[0,0,0,0,0,0,0,0,0,0],
mode='markers',
text=[ 'Battle of Marawi'
],
marker = dict(
color = '#d3d3d3',
size = 5,
opacity = 0.0,
line = dict(
color = 'black',
width = 2
)
),
showlegend = False
)
]
layout = go.Layout(
title = 'Raw Car Counts Per Scene',
xaxis = dict(
title='Date'),
yaxis = dict(
title='Number of Cars'),
# to highlight the timestamp we use shapes and create a rectangular
shapes = [
# 1st highlight during 23 May 2017 - 17 Oct 2017
{
'type': 'rect',
'layer':'below',
# x-reference is assigned to the x-values
'xref': 'x',
# y-reference is assigned to the plot paper [0,1]
'yref': 'paper',
'x0': '2017-5-23',
'y0': 0,
'x1': '2017-10-17',
'y1': 1,
'fillcolor':oic.OI_LIGHT_BLUE_HEX ,
'opacity': .5,
'line': {
'width': 0,
}
}
]
)
fig = go.Figure(data=data, layout=layout)
py.iplot(fig,filename='cufflinks/line')
The graph above shows raw car counts for the combat-area AOI. The points correspond to counts from scenes that cover at least 25% of the AOI. This threshold can be changed depending on the specific circumstances of the project.
from IPython.display import Image
from IPython.core.display import HTML
Image(url = "https://dl.dropboxusercontent.com/s/7tm5exsngecj2no/marawi_small_20140405_zoom.png?dl=0?raw=1", width = 500)
Inconsistent with the trend of higher car counts before the battle, is a 5 April 2014 count of only 100 cars. However, upon closer examination, this scene provided only partial coverage for the AOI. Also, the part of the AOI where there was coverage was in argueably a less densly populated area. These factors can be taken into account during post project analysis in order to create a more accurate portrait of what is happening on the ground.
Image(url = "https://dl.dropboxusercontent.com/s/hpt0b5993ot7tbz/marawi_small_20170410_zoom.png?dl=0raw=1", width = 500)
Cars visualized as dots from an image dated 10 April 2017. 892 cars were counted in this part of Marawi, just over a month before the battle began. This scene had total AOI coverage.
Image(url = "https://dl.dropboxusercontent.com/s/8mb39svv65jaq3i/marawi_small_20170529_zoom.png?dl=0?raw=1", width = 500)
Cars visualized as dots from an image dated 29 May 2017. 102 cars were counted in this part of Marawi a week after the battle began. The scene had total AOI coverage.
Because of the potentially randomized nature of single-collect observations, we can also aggregate detection statistics along monthly intervals. In this case, figures are summarized by the maximum and average values of monthly composites.¶
# Next will generate monthly average/max versions:
# Create a new dataframe with just the datetime, column of interest (obj_count), and city
plt_df_m = raw_df[['rounded_local_time', 'obj_count', 'city']]
# This^ will apply to all cities (when fix and add other project-ljw), don't need to repeat elsewhere
# This section is just for Marawi, we need to split them out before indexing
plt_df_m_dam = plt_df_m[plt_df_m['city']=='Marawi']
plt_df_m_dam = plt_df_m_dam.groupby(['rounded_local_time', 'city'])['obj_count'].sum().reset_index()
# Then we can set the datetime column as the index
plt_df_m_dam.index = plt_df_m_dam['rounded_local_time']
# Group the data by month, and take the max mean for each group (i.e. each month)
# Note for group assignment: M = Monthly End, MS = Monthly Start, BMS = First Business Day of Month
plt_df_m_dam_max = plt_df_m_dam.resample("MS").max()
plt_df_m_dam_mean = plt_df_m_dam.resample("MS").mean()
# Now we'll plot the monthly high and mean:
data = [
go.Bar(
x=plt_df_m_dam_max.index,
y=plt_df_m_dam_max['obj_count'],
name='Monthly High',
marker=dict(
color=oic.OI_GREEN_HEX
),
),
go.Bar(
x=plt_df_m_dam_mean.index,
y=plt_df_m_dam_mean['obj_count'],
name='Monthly Average',
marker=dict(
color= oic.OI_DARK_BLUE_HEX
),
),
go.Scatter(
x=['2017-5-23','2017-7-2','2017-09-22','2017-10-17'],
y=[1,1,1,1],
name='',
mode='markers',
text=['Battle of Marawi Begins',
'Martial Law Declared',
'AFP Retakes Masiu Bridge',
'Marawi Declared Liberated'],
marker = dict(
color = 'black',
size = 5,
opacity = 0.5,
line = dict(
color = 'black',
width = 2
)
),
showlegend = False
)
]
layout = go.Layout(
title = 'Monthly Car Counts',
xaxis = dict(
title='Date'),
yaxis = dict(
title='Number of Cars'),
# to highlight the timestamp we use shapes and create a rectangular
shapes = [
# 1st highlight during battle
{
'type': 'rect',
'layer':'below',
# x-reference is assigned to the x-values
'xref': 'x',
# y-reference is assigned to the plot paper [0,1]
'yref': 'paper',
'x0': '2017-05-23',
'y0': 0,
'x1': '2017-10-17',
'y1': 1,
'fillcolor': 'gray',
'opacity': 0.4,
'line': {
'width': 0,
}
}
]
)
fig = go.Figure(data=data, layout=layout)
py.iplot(fig,filename='cufflinks/line')
The monthly high and monthly average are helpful in a few ways. When looking at historical imagery in 2009, many monthly highs and monthly means are the same. This would indicate only one scene was taken during that month. Yet, overall we can still see a decline in car counts once the battle started. In addition, it is another way to track reconstruction and repopulation efforts post-battle, assuming higher revisit rates continue as compared to pre-battle revisit rates.
Incorporating External Data, Graphics¶
Open Source reporting can also be incorporated to contextualize findings. In this case, the aforementioned offensive correlates both temporally and spatially with declining car counter values.
While this example is derived from external open sources, we can also generate our own static or interactive graphics based on event data feeds or other interactive resources.
from IPython.display import Image
from IPython.core.display import HTML
Image(url= "https://celestialweekly.com/wp-content/uploads/2017/07/14marawi-graphic-facebookJumbo-v4-1024x536.png", width = 500)
Source: Celestial Weekly
Car Dispersion Within Combat Zone¶
In addition to the overall impact of conflict on traffic within Marawi, spatial patterns can be analyzed to derive additional conclusions or open up new lines of analytic questioning. We are able to generate an interactive webmap displaying commonly accessible background maps underneath CV-derived car points. This allows analysts to display and explore their findings in a readymade graphic format.
Prior to the start of the battle, the most densly populated car counts were found in the west/southwest part of the the conflict area AOI. Once the battle began, and civilians displaced, the general civilian POL changed. Car counts decreased overall within the combat AOI, and their dispersion argueably changed as well. This is especially evident in later stages as well as post battle.
If a major city normally has the highest car concentration in the central business district (CDB), and that city was the location of a battle, then it can be assumed that that normal civilian pattern of life would be affected. This would mean that there may not be as high of a concentration of cars in the CBD during or after the battle, when fighting and/or reconstruction is taking place.
By understanding the counts and dispersion of cars, an analyst can better identify civilian and/or military patterns of life. In this case, the removal of civilian cars in an area of the AOI that normally has the highest concentration of cars could queue the analyst to focus more on the movements of the cars that are still there. Furthermore, if a high concentration of cars returns to what was considered the normal POL car concentration location, then it could indicate a return of displaced persons.
def make_feature_group(date,color,feature_group):
dd = raw_df[raw_df['city']=='Marawi']
dd = oi_pandas.read_df(dd,geom_columns=['aoi_measured_geom'], use_shapely=False)
dd = dd[date]
cent = dd['aoi_measured_geom'].iloc[0].Centroid()
vals = []
for ix, row in dd.iterrows():
for i in range(0,row['aoi_measured_geom'].GetGeometryCount()):
y = row['aoi_measured_geom'].GetGeometryRef(i).GetX()
x = row['aoi_measured_geom'].GetGeometryRef(i).GetY()
vals.append((x,y))
dat = pd.DataFrame(vals)
dat = dat.rename(columns={0:'latitude',1:'longitude'})
for ix, row in dat.iterrows():
folium.CircleMarker(location=[row.latitude, row.longitude],
radius=5,
color=color,
weight=3).add_to(feature_group)
m = folium.Map(prefer_canvas=True,location=[7.999944, 124.301010],zoom_start=14)
feature_group = FeatureGroup(name='2017-3-14: ~2 Months Before Battle')
feature_group1 = FeatureGroup(name='2017-4-10: ~1 Month Before Battle')
feature_group2 = FeatureGroup(name='2017-5-29: 1 Week Into Battle')
feature_group3 = FeatureGroup(name='2017-9-18: ~4 Months Into Battle')
feature_group4 = FeatureGroup(name='2018-6-7: ~7 Months After End of Battle')
make_feature_group('2017-3-14',oic.OI_GRAY_HEX,feature_group)
make_feature_group('2017-4-10',oic.OI_GREEN_HEX,feature_group1)
make_feature_group('2017-5-29',oic.OI_LIGHT_BLUE_HEX,feature_group2)
make_feature_group('2017-9-18',oic.OI_PURPLE_HEX,feature_group3)
make_feature_group('2018-6-7',oic.OI_MAGENTA_HEX,feature_group4)
feature_group.add_to(m)
feature_group1.add_to(m)
feature_group2.add_to(m)
feature_group3.add_to(m)
feature_group4.add_to(m)
LayerControl().add_to(m)
m
Sample scenes in the above graphic were chosen based on the best AOI coverage and lowest cloud score. The dates reflect pre, during, and post battle timeframes. With this data, an analyst can begin to analyze the spatial distribution of cars over time and make assumptions concerning what is happening on the ground:
- 14 MAR 2017: Before battle; a normal pattern of life with car counts clustered in the CBD.
- 10 APR 2017: Before battle; another normal pattern of life scene observed, cars still clustered in CBD.
- 29 MAY 2017: One week into battle; severe decline in car counts overall, still some clustering in CBD, potential exodus of civilians.
- 18 SEP 2017: One month before end of battle; continued decline in car counts, less clustering, potential removal of most civilian vehicles.
- 07 JUN 2018: Seven months after battle; car counts beginning to increase, no CBD clustering, potential reconstruction efforts but no return of civilian POL.
Landuse/Landcover (LULC) Analysis¶
Temporal Landcover Change:¶
According to the United Nations High Commissioner for Refugees (UNHCR) 98% of the population of Marawi City was forcibly displaced.
When the battle was declared over, many of the city's buildings were either destroyed or damaged due to aerial bombing, fire, or other explosive weaponry. By running a landuse algorithm over available imagery, an analyst can begin to track LULC classification changes.
The following visuals could also help military planners have a better high-level understanding of what is happening in their area of operations. These visuals could eventually be incorporated into maps or other tactical mission planning products that are used to brief mission orders. Subsequent mission planning products could then be physically taken on a tactical mission where the status of a building or road is needed and can be easily referenced by the operator.
df_lulc = pd.read_csv("/orbital/use_cases/marawi/marawi_big_monthly/area_20170422_20170522_to_20170523_20180623_all.csv")
trace1 = go.Bar(
x=['construction'],
y=[155],
marker=dict(
color=[oic.OI_LIGHT_BLUE_HEX
]),
name='Construction Polygon Detected'
)
trace2 = go.Bar(
x=['destruction'],
y=[698],
marker=dict(
color=[oic.OI_GREEN_HEX
]),
name='Destruction Polygon Detected'
)
data = [trace1, trace2]
layout = go.Layout(
barmode = 'group',
title='Baseline: 1 Month Before Battle - Target: 13 Months After Start of Battle',
xaxis=dict(
title='Change Type',
titlefont=dict(
size=18,
color='black'
)
),
yaxis=dict(
title='Polygon Change Count',
titlefont=dict(
size=18,
color='black'
)
)
)
fig = go.Figure(data=data, layout=layout)
py.iplot(fig, filename='grouped-bar')
a = Image(url='https://dl.dropboxusercontent.com/s/4tvbfbebrwm6ko9/2016_to_2017_before.png?dl=0', width = 1000)
display(a)
Comparing two pre-battle time ranges, relatively little construction or destruction.
b = Image(url='https://dl.dropboxusercontent.com/s/sik1j245fkxij1y/april_to_june.png?dl=0', width = 1000)
display(b)
Comparing one month before to one month after the battle begins, destruction polygons are observed in heaviest area of fighting.
c = Image(url='https://dl.dropboxusercontent.com/s/2txmxcj1pgwtenc/may_to_july.png?dl=0', width = 1000)
display(c)
As fighting continues, more strategic targets are engaged. The group of destruction polygons above could indicate goverment forces closing in on a Maute terrorist stronghold. In this case, that stronghold is a mosque and houses the Masjed Mindanao Islamic Centre.
"The mosque, which houses the Masjed Mindanao Islamic Centre, had been used by militants from the Maute terrorist group as a shelter, stockade, sniper's nest and holding area for their hostages. Satellite images show that the area around it has been reduced to rubble. The army had avoided bombing it out of respect for Muslims." - The Strait Times (https://www.straitstimes.com/asia/govt-troops-retake-grand-mosque-in-marawi)
Concerning the construction polygons of roads, this could be due to an image in the before time range with a road covered with rubble, but then that rubble being cleared in the after time range. It would make sense that rubble would be routinely created by bombings, and then removed to some degree in order to facilitate mounted patrols or other vehicle-necessary military operations.
j = Image(url='https://dl.dropboxusercontent.com/s/ja0jh65ghb66zu1/2070711_01_02_mosque1.png?dl=0', width = 800)
display(j)
Above shows a close up of changes near the Masjed Mindanao Islamic Centre. These changes were detected using a LULC algorithm trained on medium resolution Planet imagery.
i = Image(url='https://dl.dropboxusercontent.com/s/kuv31o18nlpoaz1/20170712_DG_mosque1.png?dl=0', width = 1000)
display(i)
Above is a high resolution Digital Globe image taken a day after the landuse algorithm detected change in Planet imagery. Even within commercial imagery sources, a variation of "tip and cue" could be used.
d = Image(url='https://dl.dropboxusercontent.com/s/2qpiz5u77tu5kba/june_to_august.png?dl=0', width = 1000)
display(d)
e = Image(url='https://dl.dropboxusercontent.com/s/m8jmqy3e9xo83nd/july_to_september.png?dl=0', width = 1000)
display(e)
f = Image(url='https://dl.dropboxusercontent.com/s/n2wk3c9fqwj1o8e/aug_to_oct.png?dl=0', width = 1000)
display(f)
g = Image(url='https://dl.dropboxusercontent.com/s/314r7gkaech4rps/september_to_november.png?dl=0', width = 1000)
display(g)
k = Image(url='https://dl.dropboxusercontent.com/s/86vgcglzs5qim00/pl_constr_helipad.png?dl=0', width = 300)
display(k)
Construction can also be examined in the AOI during the battle. Above, construction polygons were added as month to month analysis was observed.
l = Image(url='https://dl.dropboxusercontent.com/s/pvb4nx7a6bvlsbk/dg_constr_helipad.png?dl=0', width = 1000)
display(l)
By examining the month to month construction polygons, a developing pattern would urge an analyst for further investigation. When looking at Digital Globe high resolution imagery, we find a helipad. Together, the helipad and construction of new buildings could indicate a field hospital, or other AFP military installation. By identifying other such installations, an order of battle can be pieced together. Detecting this trend without the help of the LULC algorithm would require increased analyst time and resources.
h = Image(url='https://dl.dropboxusercontent.com/s/nhfnq3rlz1f1kxh/pre_1month_to_oct.png?dl=0', width = 1000)
display(h)
df_lulcdes = pd.read_csv("/orbital/use_cases/marawi/marawi_big_monthly/area_destruction_20170422_20170522_to_20170523_20180623_all.csv")
data = [
go.Scatter(
x=df_lulcdes['change_ts'],
y=df_lulcdes['area_m2_run_tot'],
name='Marawi',
mode = 'lines+markers',
line = dict(color = (oic.OI_DARK_BLUE_HEX),
width = 2),
marker = dict(color = (oic.OI_DARK_BLUE_HEX),
size = 5)
),
go.Scatter(
x=['2017-08-10'],
y=[0,0,0,0,0,0,0,0,0,0],
mode='markers',
text=[ 'Battle of Marawi'
],
marker = dict(
color = '#d3d3d3',
size = 5,
opacity = 0.0,
line = dict(
color = 'black',
width = 2
)
),
showlegend = False
)
]
layout = go.Layout(
title = 'Cumulative Destruction Area',
xaxis = dict(
title='Date'),
yaxis = dict(
title='Square Meters'),
# to highlight the timestamp we use shapes and create a rectangular
shapes = [
# 1st highlight during 23 May 2017 - 17 Oct 2017
{
'type': 'rect',
'layer':'below',
# x-reference is assigned to the x-values
'xref': 'x',
# y-reference is assigned to the plot paper [0,1]
'yref': 'paper',
'x0': '2017-5-23',
'y0': 0,
'x1': '2017-10-17',
'y1': 1,
'fillcolor':oic.OI_LIGHT_BLUE_HEX ,
'opacity': .5,
'line': {
'width': 0,
}
}
]
)
fig = go.Figure(data=data, layout=layout)
py.iplot(fig,filename='cufflinks/line')
df_lulcdes = pd.read_csv("/orbital/use_cases/marawi/marawi_big_monthly/area_destruction_20170422_20170522_to_20170523_20180623_all.csv")
df_lulcdes = pd.read_csv("/orbital/use_cases/marawi/marawi_big_monthly/area_construction_20170422_20170522_to_20170523_20180623_all.csv")
data = [
go.Scatter(
x=df_lulcdes['change_ts'],
y=df_lulcdes['area_m2_run_tot'],
name='Marawi',
mode = 'lines+markers',
line = dict(color = (oic.OI_DARK_BLUE_HEX),
width = 2),
marker = dict(color = (oic.OI_DARK_BLUE_HEX),
size = 5)
),
go.Scatter(
x=['2017-08-10'],
y=[0,0,0,0,0,0,0,0,0,0],
mode='markers',
text=[ 'Battle of Marawi'
],
marker = dict(
color = '#d3d3d3',
size = 5,
opacity = 0.0,
line = dict(
color = 'black',
width = 2
)
),
showlegend = False
)
]
layout = go.Layout(
title = 'Cumulative Construction Area',
xaxis = dict(
title='Date'),
yaxis = dict(
title='Square Meters'),
# to highlight the timestamp we use shapes and create a rectangular
shapes = [
# 1st highlight during 23 May 2017 - 17 Oct 2017
{
'type': 'rect',
'layer':'below',
# x-reference is assigned to the x-values
'xref': 'x',
# y-reference is assigned to the plot paper [0,1]
'yref': 'paper',
'x0': '2017-5-23',
'y0': 0,
'x1': '2017-10-17',
'y1': 1,
'fillcolor':oic.OI_LIGHT_BLUE_HEX ,
'opacity': .5,
'line': {
'width': 0,
}
}
]
)
fig = go.Figure(data=data, layout=layout)
py.iplot(fig,filename='cufflinks/line')
def make_feature_group(date,color,feature_group):
dd = raw_df[raw_df['city']=='Marawi']
dd = oi_pandas.read_df(dd,geom_columns=['aoi_measured_geom'], use_shapely=False)
dd = dd[date]
cent = dd['aoi_measured_geom'].iloc[0].Centroid()
vals = []
for ix, row in dd.iterrows():
for i in range(0,row['aoi_measured_geom'].GetGeometryCount()):
y = row['aoi_measured_geom'].GetGeometryRef(i).GetX()
x = row['aoi_measured_geom'].GetGeometryRef(i).GetY()
vals.append((x,y))
dat = pd.DataFrame(vals)
dat = dat.rename(columns={0:'latitude',1:'longitude'})
for ix, row in dat.iterrows():
folium.CircleMarker(location=[row.latitude, row.longitude],
radius=5,
color=color,
weight=3).add_to(feature_group)
Main LULC Takeaways:¶
- Automatically derived scene to scene changes can empower an analyst to look more closely at specific landcover classes as new imagery is acquired.
- When comparing landuse output at scale before and after the battle, or month to month, the building and road classes can identify areas requiring further analysis.
- This AOI can be examined long term in order to evaluate the speed of which reconstruction takes place.
- By understanding the capabilities of automated landuse mapping, an analyst can incorporate relevant information into his/her analysis products or models.
Geolocation Data¶
As part of an ongoing partnership to acquire global geolocation data, Orbital Insight was provided a 1-month global sample. Although the sample timeframe is not ideal in regards to the time of the battle, it does showcase what possibilities there are when incorporating geolocation data:
- Identify a high level POL change within a city due to an acute deacrease of counts.
- Identify movements of distinct IDs from the confict area to another. If a high percentage of specific IDs from a conflict area are found together in another area, this could indicate a displaced persons camp.
- Understand if bombed out regions of an area are associated with "cold zones" where little to none pings are present.
- Monitor reconstruction and civilian repopulation post battle.
#! pip install s3fs
import os
from IPython.display import Image
from IPython.display import IFrame
#! pip install s3fs
import s3fs
geo_raw = pd.read_csv('s3://orbital-qubole-output-store/logan/marawi/data_all/marawi_big/part-00000-94e70d24-092a-4dda-b81f-175014f70a06-c000.csv')
geo_raw.set_index('date', inplace=True)
n = geo_raw['ad_id'].count()
u = geo_raw['ad_id'].nunique()
#print('Total number of pings',n)
#print ('Total Number of unique devices',u)
geo_raw1 = pd.read_csv('s3://orbital-qubole-output-store/logan/marawi/data_all/iligan_detail/part-00000-41d52146-e7d7-4ec2-9b7f-06fece1e4a0a-c000.csv')
geo_raw1.set_index('date', inplace=True)
i = geo_raw1['ad_id'].count()
t = geo_raw1['ad_id'].nunique()
#print('Total number of pings',i)
#print('Total number of unique devices',t)
# Hexbin Plotter
# Path to Input CSV
inputCSV = "s3://orbital-qubole-output-store/logan/marawi/data_all/marawi_big/part-00000-94e70d24-092a-4dda-b81f-175014f70a06-c000.csv"
# Area WKT
areaWKT = "POLYGON ((124.2693312266717385 8.0242587437599013, 124.3328492909273564 8.0239620636372138, 124.3320503215704207 7.9828203268138367, 124.2434658377126198 7.9832159401960041, 124.2439651935557094 8.0238631702215173, 124.2693312266717385 8.0242587437599013))"
# hexbinHeight (Optional)
hexbinHeight = 0.001
# mapStyle
mapStyle = "OpenStreetMap"
# Read CSV
csv_df = pd.read_csv(inputCSV, encoding = 'utf-8')
# Make Hex Bins Using the Script
hexBinMaker = HexBinMaker(csv_df, areaWKT, hexbinHeight, verbose = False)
hexbinCounts = hexBinMaker.makeHexBins()
# Plot Map
m = hexBinMaker.plotHexBins(hexbinCounts, areaWKT, vmin = 1, vmax = 300)
m.save("marawi_hexbin_example.html")
m
Between June 2018 and July 2018, geolocation pings for Marawi:
- Total number of pings: 117,336
- Total Number of unique devices: 4,099
Note that despite the battle being over for about 8 months, geolocation pings in the southeast part of the city are still reduced.
With more expansive geolocation data, monitoring the return of civilians and/or estimating reconstruction efforts can be analyzed.
# Hexbin Plotter
# Path to Input CSV
inputCSV = 's3://orbital-qubole-output-store/logan/marawi/data_all/iligan_detail/part-00000-41d52146-e7d7-4ec2-9b7f-06fece1e4a0a-c000.csv'
# Area WKT
areaWKT = 'POLYGON ((124.2186874635650042 8.2094699920194500, 124.2187302757100014 8.2058258591004396, 124.2502828344649970 8.2063767185415504, 124.2583482350650002 8.2122023547542096, 124.2628387421169975 8.2207993748828301, 124.2673971779610014 8.2269987151005495, 124.2812660454069942 8.2465329811295902, 124.2795417958029986 8.2604208378065902, 124.2747596891450002 8.2581824898584593, 124.2692667287130064 8.2576708656947897, 124.2639030144780037 8.2558162225211600, 124.2589270385749955 8.2543452924582699, 124.2568591005580032 8.2515313239738592, 124.2558251315499973 8.2499324693269500, 124.2513029876180042 8.2499324692918794, 124.2498812802249972 8.2509557370306794, 124.2481886606140051 8.2537782685039893, 124.2469943108870041 8.2539095994244498, 124.2466957234559999 8.2537126030142591, 124.2440558109170041 8.2517323673459906, 124.2434414263070011 8.2511919012833896, 124.2429294391329933 8.2495367193491695, 124.2428099754500010 8.2488611329082495, 124.2429806378229955 8.2479490893567498, 124.2425539818350018 8.2459898776109100, 124.2418030673489966 8.2445035725591502, 124.2412740139269971 8.2436084088645991, 124.2409122781230053 8.2429183923342304, 124.2408923313020068 8.2424989013711105, 124.2407876105559978 8.2420448635395207, 124.2405981158519950 8.2409295076120106, 124.2403587541290051 8.2405001440714791, 124.2393602383189943 8.2392311344245002, 124.2387618340419948 8.2384957842308495, 124.2382979176300069 8.2373020088691593, 124.2381433299469933 8.2370108080699804, 124.2375548990480070 8.2366110518796098, 124.2366692934070045 8.2357193553506498, 124.2347672888650010 8.2346090821983307, 124.2336703703419971 8.2351275857511297, 124.2322410196020002 8.2323196223221906, 124.2330667401219984 8.2318292948193594, 124.2332603574210026 8.2306016664557795, 124.2334255015179991 8.2288581316710907, 124.2330667401650004 8.2286552358288692, 124.2330667401210036 8.2282564007905297, 124.2344612971400011 8.2281360971416895, 124.2344385186300002 8.2276288564350004, 124.2336657279250005 8.2267815528425707, 124.2329254267459930 8.2262010419901799, 124.2320450264700042 8.2229145262764494, 124.2326882926199971 8.2219064914073101, 124.2313481548220011 8.2186701516037495, 124.2293647508739980 8.2168662786079700, 124.2274885579650032 8.2161765603038592, 124.2242722272579982 8.2130993409089701, 124.2213239240840039 8.2103404342966702, 124.2186874635650042 8.2094699920194500))'
# hexbinHeight (Optional)
hexbinHeight = 0.001
# mapStyle
mapStyle = "OpenStreetMap"
# Read CSV
csv_df = pd.read_csv(inputCSV, encoding = 'utf-8')
# Make Hex Bins Using the Script
hexBinMaker = HexBinMaker(csv_df, areaWKT, hexbinHeight, verbose = False)
hexbinCounts = hexBinMaker.makeHexBins()
# Plot Map
m = hexBinMaker.plotHexBins(hexbinCounts, areaWKT, vmin = 1, vmax = 300)
m.save("iligan_hexbin_example.html")
m
Between June 2018 and July 2018, geolocation pings for Iligan:
- Total number of pings: 461,934
- Total Number of unique devices: 18,624
Unlike Marawi, Iligan shows what can be perceived as a normal city's geolocation distribution.
